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Abstract 



The dynamics of a self-gravitating cold Fermi gas is described using the analogy 

with an interacting self-gravitating Bose condensate having the same Thomas-Fermi 

limit. The dissipationless formation of a heavy neutrino star through gravitational 

collapse and ejection of matter is demonstrated numerically. Such neutrino stars 

ON I offer an alternative to black holes for the supermassive compact dark objects at the 

^ ■ centers of galaxies. 

CN ■ 
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Supermassive neutrino stars, in which self-gravity is balanced by the degener- 
acy pressure of the cold fermions, have been a subject of speculation for more 
than three decades Originally, these objects were proposed as models for 
dark matter in galactic halos and clusters of galaxies, with neutrino masses 
in the ~ eV range. More recently, however, degenerate superstars composed 
of weakly interacting fermions in the ~ 10 keV range, have been suggested as 
an alternative to the supermassive black holes that are purported to exist at 
the centers of galaxies In fact, it has been shown [|J that such degen- 

erate fermion stars could explain the whole range of supermassive compact 
dark objects which have been observed so far, with masses ranging from 10 6 
to 3 x 10 9 M Q , merely assuming that a weakly interacting quasi-stable fermion 
of mass rrif ~ 15 keV exists in nature. 

As an example, the most massive compact dark object ever observed, is located 
at the center of M87, with a mass M ~ 3.2 x 1O 9 M |7[. Interpreting this object 
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as a relativistic fermion star at the Oppenheimer-Volkoff || limit, yields a 
fermion mass of rrif ~ 15 keV and a fermion star radius of R = 4.45i?s — 1-5 
light-days |||4|], where R$ is the Schwarzschild radius. In this case there is 
little difference between the fermion star and black hole scenarios, because the 
radius of the last stable orbit around a Schwarzschild black hole is R = 3Rs 
anyway. 



Extrapolating this down to the compact dark object at the center of our galaxy 
||, which, having a mass M ~ 2.6 x 1O 6 M , is at the lower limit of the mass 
range of the observed compact dark objects, we obtain, using the same fermion 
mass, a radius R ~ 20 light-days ~ 5 x 10 4 i?s 0. As the potential inside such 
a nonrelativistic fermion star is shallow, the spectrum of radiation emitted by 
accreting baryonic matter is cut off for frequencies larger than 10 13 Hz ||||, 
as is actually observed in the spectrum of the enigmatic radio source Sgr A* 



at the galactic center jnj. A fermion star with radius R < 20 light-days and 
mass M ~ 2.6 x 1O 6 M is also consistent || with the observed motion of stars 
within a projected distance of 6 to 30 light-days from Sgr A* M. 



Of course, it is well-known that 15 keV lies squarely in the cosmologically 
forbidden mass range for stable active neutrinos v ||11|| . The existence of such a 



massive active neutrino is also disfavoured by the Super-Kamiokande data [12 



However, as shown by Shi and Fuller ]T3| for an initial lepton asymmetry of 
~ 10~ 3 , a sterile neutrino v s of mass m s ~ 10 keV may be resonantly produced 
in the early universe with near closure density, i.e., Q ~ 1. The resulting 
energy spectrum is not thermal but rather cut off so that it approximates a 
degenerate Fermi gas. Sterile neutrinos in this mass range are also constrained 
by astrophysical bounds on the radiative decay v s — > z/7 |Tj|. 
allowed parameter space includes m s ~ 15 keV, contributing Q 



However, the 
, ~ 0.3 to the 



critical density, as favoured by the BOOMERANG data [15 



The statics of degenerate fermion stars is well understood, being the Op- 
penheimer-Volkoff equation in the relativistic case pi, or the Lane-Emden 
equation with polytropic index n = 3/2 in the nonrelativistic limit [plj. Al- 
ternatively, one may understand these as the Thomas-Fermi theory applied to 
self-gravitating fermion systems. The extension of the Thomas-Fermi theory 
to finite temperature |T7],|23j has been used to show that, at a certain critical 
temperature, weakly interacting massive fermionic matter undergoes a first- 
order gravitational phase transition from a diffuse to a clustered state, i.e. 
a nearly degenerate fermion star. However, such studies do not bear on the 
crucial dynamical question of whether the fermion star can form through grav- 
itational collapse of density fluctuations in an orthodox cosmological setting. 
Indeed, since collisional damping is negligible, one would expect that only a 
virialized cloud results [ITT] . 



iV-body simulations of the collisionless Boltzmann or Vlasov equation evidence 
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a rather different picture: the collapse is followed by a series of bounces, with 



matter expelled at each bounce, leaving behind a condensed object Hl9fl . By 
Liouville's theorem, the Vlasov equation describes an incompressible fluid in 
phase-space so that it respects a form of the exclusion principle. Hence, these 
iV-body simulations are effectively fermion simulations. What transpires is 
that gravity, being attractive, self-organises the phase-space fluid into a high- 
density /momentum core at the expense of other low-density/momentum re- 
gions as seen in the evolution of the spherical Vlasov equation . 

Much the same behaviour is observed in the formation of mini-boson stars 



through so-called gravitational cooling ||21|| . Such a boson star is stable by 
balancing the uncertainty and gravitational pressures. A similar mechanism 
works in the presence of a quartic self-interaction p2|, which dominates over 



uncertainty pressure, resulting in an equilibrium radius R 3> 1/m, where m 
is the boson mass E3]. Hence, we have a universal description of the physics 



underlying the formation process: once the collapse proceeds far enough, un- 
certainty, interaction or degeneracy pressure results in a bounce, in which the 
outgoing shock wave carries away the binding energy. The virial counter ar- 
gument mentioned above is bypassed, because the ejected matter invalidates 
its assumption that there is no flow through the boundary. 

The purpose of this letter is to verify this picture for the formation of a 
fermion star from a cold gravitationally unstable configuration. The dynam- 
ical Thomas-Fermi theory, developed long ago by Bloch for the electron gas 
24|, amounts to Euler's equations for irrotational flow, subject to an equa- 



tion of state P = -P(p). The problem is that, due to the presence of shocks 
and instabilities, a naive integration of the Bloch equations is precluded in 
the gravitational case. The usual remedy is to introduce some small numerical 
viscosity. However, it seems imprudent to draw conclusions based on introduc- 
ing an ad-hoc dissipation into what is fundamentally a dissipationless system. 
Here we take another, literally conservative approach. 

We base our approach on the equivalence of a degenerate non-interacting 
fermion star with a certain type of self-interacting cold boson star. We shall 
demonstrate this equivalence using scaling arguments which are similar in 
spirit to the analysis of self-interacting boson stars by Colpi, Shapiro and 
Wassermann Consider a complex scalar field ip with a repulsive La- 
grangean interaction term U(\ip\ 2 ). For simplicity, we use the Newtonian ap- 
proximation, but our analysis may be easily extended in a general relativistic 
context. In the Newtonian limit, a self-interacting boson star is governed by 
the Gross-Pitaevskii-like equations 



dip 



A 1 dP 



1>, (1) 



2m m d\ip\ r " 

Ay? = 4vrGp, (2) 
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p = m 2 \ip\ 2 . 



(3) 



The pressure tensor is given by 



_ dip dip* 

Re ^~~ a % 
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A| 



(4) 



Following Colpi et al. [23|, we introduce a dimensionless parameter A and 
define 



R* = — = A— ^; M* = i?*mp. 



m 



(5) 



as the length and mass scales, respectively. Here mpj = 1/y/G denotes the 
Planck mass and the arbitrary constant A will be fixed later using physical 
arguments. The substitution 



iP 



m 



4ttA 



yields the coupled dimensionless equations 
A 



% d^> 
A ~dt ~ 



2A 2 



+ <p + V(\^\ 



(6) 



(7) 



where we have introduced the dimensionless potential 

U l) ~ m 4 d|^| 2 ' 
The pressure tensor M), written in terms of dimensionless variables, reads 



(9) 



Pi, 



m 1 
i^A 2 A 2 



Re 



d^ d^* 

dxj dxn 
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1^1 



dU 
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(10) 



A static, spherically symmetric solution, usually referred to as boson star, is 
obtained by the ansatz 



$(r), 



;n) 



where $(r) is a real function. It is clear from the definition fl5|) that A ^> 1, as 
long as m/ropi <C A which is a reasonable assumption. Hence, for sufficiently 
large A, Eq. (^) with ([TTD degenerates to 



m 



0. 



(12) 



This equation is exact in the limit A — > oo, which is just the Thomas- Fermi 
limit [25]. In this limit the derivative terms in Eq. fliTf ) vanish so that the 
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pressure tensor becomes diagonal, with all the components equal to P = Pa. 
We obtain an effective equation of state given by 



m 



■i 



»=4^*' (13) 
P = pV -U. (14) 



It may be easily shown that the last equation combined with ([121 ) yields the 
equation of hydrostatic equilibrium 

dP 

— = -dip. (15) 

For a particular U(\ip\ 2 ), the equation of state in the form P = P(p) is obtained 
by eliminating $ 2 from (|1~3|) and (0). On the other hand, if the equation of 
state, e.g., of the polytropic type, is given, we can determine the potential U 



and its derivative V by integrating ( [T5|) . For a general polytropic equation of 
state 

P(p) = Kp 1+1/n , (16) 
Eq. (H) yields 

V = (n + l)Kp 1/n . (17) 

This, together with Eq. ([03), gives the potential in terms of $ 2 . As we still 
have the freedom to choose the parameter A conveniently, we fix A so that the 
potential V in (0) takes a simple form 

V=\^\ 2/n . (18) 

The polytropic equation of state with n = 3/2, together with the equation of 
hydrostatic equilibrium ([TBI) and Poisson's equation @, describes a degener- 
ate fermion star. Hence, we have demonstrated that a degenerate fermion star 
is equivalent to a self- interacting boson star in the limit A — > oo. Moreover, it 
may be shown numerically that, even for moderate values of A of the order of 
10 to 20, the static solutions are almost degenerate and are quite well approx- 
imated by the static solution for an infinite A. This has been demonstrated for 
a quartic self-interaction p3| (note that our A corresponds to a/A in reference 



Next we proceed to solving Eqs. (|7|), (|8|) and ([IB] ) numerically for large, but fi- 
nite value of A. For weakly interacting degenerate fermions with the polytropic 
index n = 3/2, the length and mass scales are given by 

/9tt 2 \ 1/4 mpi /keV\ 2 [2 , s 

R * = ^ = - 2325 \ ~ ^ ' 19 
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M* = 1.489 (—} J^- x 10 12 M , (20) 
\ m f J \ 9f 

where gf is the spin degeneracy and rrif the mass of the fermion. As M* is 
of the order of the Oppenheimer-Volkoff limit, the validity of the Newtonian 
approximation in the static case requires 

M = (An)- 1 J d 3 r |^| 2 < 1. (21) 

By construction, once the overall scale is specified, a large but finite A allows 
us to simulate the fermionic problem as a bosonic one, as long as their Thomas- 
Fermi limits coincide, while providing an explicitly energy conserving way of 
controlling the shocks and instabilities. The basic regulating mechanism is the 
kinetic part of (0) which penalises density spikes. Of course, A must be so 



large that this term does not change the static scaling relationship jig 

MR^)=C n , (22) 



arising from the polytropic equation of state. Our criterion is that the ratio 
of kinetic and pressure contributions to the static energy functional should be 
small. In particular for a Gaussian \1/ = aexp [— (r / f3) 2 ] we find 

1 Jd 3 r|Vfr| 2 / IN 3 / 2 3 
2A 2 / d 3 r|^| 2 + 2 /™ \ n) 2A 2 f3 2 a 2 / n ' 1 ' 

For n = 3/2 this is independent of the size of /3 for a given mass, yielding 
the weak condition M ^> 0.91A" 3 . When the mass violates this inequality the 
self-interaction V becomes irrelevant and there is a crossover to mini-boson 
star behaviour MR = const. 

In Fig. |l] we display the evolution of |r\l/| 2 for the spherical collapse of a mass 
M = 0.008, initially in the form of the above Gaussian with f3 = 100 = A. 
A fairly coarse sampling is used to make this contour plot, as this allows 
one to see the basic features of matter ejection without the complicated fine 
detail hiding it. Eqs. (|?|) and (Q) are first written in terms of rj = r^f, and 
then solved assuming spherical symmetry using the method of lines for 



r = to 1000. The region is divided into 4000 intervals, and a fourth order 
polynomial is fitted in each interval, with the function continuous up to its 
second derivatives across each interval. The time step is adjusted dynamically 
to achieve the desired accuracy of 10~ 6 . The boundary conditions at r = 
and r = 1000 are reflective. To prevent matter from being artificially reflected 
by the boundary at r = 1000, we have introduced an r-dependent imaginary 



part to the potential, or 'sponge' p2| , from r = 700 to 1000, which removes 
the ejected fermionic matter. The mass in the region remains constant to 
within the numerical integration error until t = 5 x 10 4 , which is a quarter 
of the simulation period, when the first ejected matter reaches the 'sponge'. 
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The expected features of bounce and ejection, leaving a condensed core, are 
evident. We obtain identical results using the Crank-Nicholson method. 



A pertinent question is whether the Thomas-Fermi approximation is satisfied 
during the ejection process. Similar to WKB, the Thomas-Fermi limit is ob- 
tained when the scale of variation of density is small compared to that of the 
potential V, i.e., 



|Vp| 2 


|V|^|| 2 


p 2 8 dU/d\ip\ 2 


2A|^|2+2/n 



r = jnrj^hnr, = „ J.^/n « L ( 24 ) 



In Fig. ^| we exhibit \r^\ 2 and £ versus r at time 1.22 x 10 4 , where ejecta first 
develop. Evidently, the Thomas-Fermi condition, is well satisfied, except near 
sharp density minima. This is expected since in a pure Thomas-Fermi system 
evolving under Bloch's equations, the sharp density gradient would evolve into 
a shock, which violates the Thomas-Fermi condition. 

In Fig. ^| we show |r\l/| 2 on the time slice t = 2 x 10 5 of our A = 100 simulation. 
The core mass M = 0.0057 and radius R = 28 are roughly commensurate with 
(22) and [16 C3/2 = 132.3843. Of course, at this point there is still evolution 



and oscillation of the core, similar to the mini-boson star case [21], with R 
varying between about 25 and 31. In general one can expect a persistent 
seismology with a period T of the order of the radius divided by the average 
speed of sound, T ~ 51/M, which agrees with the simulations. The seismology 
includes higher-order modes as well. 

In summary, using a bosonic representation of the dynamical Thomas-Fermi 
theory for a self-gravitating gas, we have shown that nonrelativistic, degener- 
ate and weakly interacting fermionic matter will form supermassive fermion 
stars through gravitational collapse accompanied by ejection. For a fermion 
mass of rrif ~ 15 keV, and a final mass M ~ 2.6 x 1O 6 M such a superstar 
is consistent with the observations of the compact dark object at the cen- 
ter of our galaxy. A similar demonstration for the formation of such a star 
near the Oppenheimer-Volkoff limit, and the question of cosmology with de- 
generate dark matter, requires a general relativistic extension which is under 
development and will be reported elsewhere. 



Acknowledgements 



The authors wish to thank Duncan Elliott (f) for many useful discussions 
regarding the simulations. This research is in part supported by the Foun- 
dation of Fundamental Research (FFR) grant number PHY99-01241 and the 
Research Committee of the University of Cape Town. The work of N.B. is 



7 



supported in part by the Ministry of Science and Technology of the Republic 
of Croatia under Contract No. 00980102. 



References 



[1] M.A. Markov, Phys. Lett. 10, 122 (1964); R. Cowsik and J. McClelland, 
Astrophys. J. 180, 7 (1973); R. Ruffini, Lett. Nuovo Cim. 29, 161 (1980). 



[2] 
[3] 

[4] 
[5] 

[e: 

[7] 
[8] 
[9 

[10 
[11 

[12 
[13 

[14 
[15 
[16 

[17 
[18 



R.D. Viollier, D. Trautmann and G.B. Tupper, Phys. Lett. B306, 79 (1993). 

N. Bilic, D. Tsiklauri and R.D. Viollier, Prog. Part. Nucl. Phys. 40, 17 (1998); 
N. Bilic and R.D. Viollier, Nucl. Phys. (Proc. Suppl.) B66, 256 (1998). 

N. Bilic, F. Munyaneza and R.D. Viollier, Phys. Rev. D59, 024003 (1999). 

D. Tsiklauri and R.D. Viollier, Astropart. Phys. 12, 199 (1999). 

F. Munyaneza, D. Tsiklauri and R.D. Viollier, Astrophys. J. 509, L105 (1998); 
ibid. 526, 744 (1999); F. Munyaneza and R.D. Viollier, |astro-ph/0103466| . 

F. Macchetto et al., Astrophys. J. 489, 579 (1997). 

J.R. Oppenheimer and G.M. Volkoff, Phys. Rev. 55, 374 (1939). 

A. Eckart and R. Genzel, Mon. N.R. Astron. Soc. 284, 576 (1997); A.M. Ghez, 

B. L. Klein, M. Morris and E.E. Becklin, Astrophys. J. 509, 678 (1998). 

R. Mahadevan, Nature 394, 651 (1998). 

E. W. Kolb and M.S. Turner, The Early Universe (Addison- Wesley, San 
Francisco, 1989). 

S. Fukuda et al., Phys. Rev. Lett. 85, 3999 (2000). 

X. Shi and G.M. Fuller, Phys. Rev. Lett. 82, 2832 (1999); K. Abazajian, G.M. 
Fuller, and M. Patel, |astro-ph/0101524. 



M. Drees and D. Wright, [hep-ph/0006274 . 

P. de Bernardis et al., Nature 404, 955 (2000). 

S. Chandrasekhar, Stellar Structure (University of Chicago Press, Chicago, 
1939). 

N. Bilic and R.D. Viollier, Phys. Lett. B408, 75 (1997); Gen. Rel. Grav. 31, 
1105 (1999); Eur. Phys. J. Cll, 173 (1999). 

W. Thirring, Z. Physik 235, 339 (1970); P. Hertel, H. Narnhofer and W. 
Thirring, Comm. Math. Phys. 28, 159 (1972); J. Messer, J. Math. Phys. 22, 
2910 (1981). 



8 



[19] J. Binney and S. Tremaine, Galactic Dynamics (Princeton University Press, 
Princeton, New Jersey, 1987), and references cited therein. 

[20] R.N. Henrikson and L.M. Widrow, Phys. Rev. Lett. 78, 3426 (1997). 

[21] E. Seidel and W.-M. Suen, Phys. Rev. Lett. 72, 2516 (1994). 

[22] J. Balakrishna, E. Seidel and W.-M. Suen, Phys. Rev. D58, 104004 (1998). 

[23] M. Colpi, S.L. Shapiro and I. Wasserman, Phys. Rev. Lett. 57, 2485 (1986). 

[24] F. Bloch, Z. Physik 81, 363 (1933). 

[25] In the context of Bose condensates see: A.S. Parkins and D.F. Walls Phys. Rep. 
303, 1 (1998). 

[26] W. Schiesser, The Numerical Method of Lines (Academic Press, San Diego, 
1991). 



9 



o 




4xl0 4 8xl0 4 1.2xl0 5 1.6xl0 5 2xl0 5 2.4xl0 5 2.8xl0 5 

Time 



Fig. 1. Combined contour-density plot for the evolution of |r\P| 2 from the initial 
configuration described in the text. Green contour lines denote levels from 10 -5 to 
10~ 4 while red lines denote levels below 1CP 5 . Gravitational collapse is followed by 
ejection of excess matter leaving a fermion star at the center. 
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